clear all; close all;

specieslist = cell(10,1);

specieslist{1} = 'Abalone';
specieslist{2} = 'Arbacia';
specieslist{3} = 'Bull';
specieslist{4} = 'Ciona';
specieslist{5} = 'Human-high';
specieslist{6} = 'Human-low';
specieslist{7} = 'Lpictus';
specieslist{8} = 'Lvariegatus';
specieslist{9} = 'Spurpuratus';
specieslist{10} = 'Zebrafish';

colors{1} = [0 1 1];
colors{2} = [1 0 0]; colors{2} = [0.33 0 0];
colors{3} = [0 0 1]; colors{3} = [0 0.5 1];
colors{4} = [1 0 1];
colors{5} = [0 0 1]; colors{5} = [0.5 0 1];
colors{6} = [0 0 1]; 
colors{7} = [1 0 0]; colors{7} = [0.67 0 0];
colors{8} = [1 0 0]; colors{8} = [1 0 0];
colors{9} = [1 0 0]; colors{9} = [1 0.5 0];
colors{10} = [0 1 0];

actualSpecies{1} = 'H. rufescens';
actualSpecies{2} = 'A. punctulata';
actualSpecies{3} = 'B. taurus';
actualSpecies{4} = 'C. intestinalis';
actualSpecies{5} = 'H. sapiens (viscous)';
actualSpecies{6} = 'H. sapiens';
actualSpecies{7} = 'L. pictus';
actualSpecies{8} = 'L. variegatus';
actualSpecies{9} = 'S. purpuratus';
actualSpecies{10} = 'D. rerio';

symb = 'oooopvs^vo';

datadate = '200617B';


%% Figure S4 plot
plotInds = [3 6 5];

ix = 0;
for k=plotInds
    flip1 = strcmp(specieslist{k}, 'Bull') || strcmp(specieslist{k},'Human-high');%strcmp(specieslist{k}, 'Bull');% || strcmp(specieslist{k},'Zebrafish')
    flip2 = 0;
    switch12 = strcmp(specieslist{k},'Human-low');%strcmp(specieslist{k},'Human-low');
    
    if flip1
       f1 = -1; 
    else
       f1 = 1;
    end
    if flip2
       f2 = -1; 
    else
       f2 = 1;
    end
    
    ix = ix+1;
    load([specieslist{k} '_' datadate '.mat']);
    
    figure(30);
    set(gcf,'pos',[900 100 320 900]);
    subplot(3,1,1);
    dm = size(U_temp{1},1);
    idx = 1:dm:(dm*(length(U_temp)-1)+1);
    Mode1 = cell2mat(U_temp);
    if switch12
        Mode2 = Mode1(:,idx);
        Mode1 = Mode1(:,idx+1);
    else
        Mode2 = Mode1(:,idx+1);
        Mode1 = Mode1(:,idx);
    end
    std1 = std(abs(Mode1),[],2);
    sOverL = (0:(length(std1)-1))/(length(std1)-1);
    med1 = median(Mode1,2);
    shadedErrorBar(sOverL,med1*f1,std1,'lineProps',{'-','Color',colors{k}});%{'r-o','markerfacecolor','r'});
    hold on;
    %plot(2,med1(2),symb(k),'Color',colors{k},'MarkerFaceColor',colors{k}); 
    ylabel('Curvature Mode 1, u_1 (a.u.)')
    xlim([0 1]); ylim([-0.9 0.9]); set(gca,'YTick',[-0.8 -0.6 -0.4 -0.2 0 0.2 0.4 0.6 0.8]);
    
    subplot(3,1,2);
    med2 = median(Mode2,2);
    std2 = std(abs(Mode2),[],2);
    shadedErrorBar(sOverL,med2*f2,std2,'lineProps',{'-','Color',colors{k}});%{'r-o','markerfacecolor','r'});
    hold on;
    %plot(2,med2(2),symb(k),'Color',colors{k},'MarkerFaceColor',colors{k}); 
    ylabel('Curvature Mode 2, u_2 (a.u.)')
    xlabel('Flagellar Arclength, s/L')
    xlim([0 1]); ylim([-0.9 0.9]); set(gca,'YTick',[-0.8 -0.6 -0.4 -0.2 0 0.2 0.4 0.6 0.8]);
    
    for i = (kfolds(1)*kfolds(2)):-1:1
        sv = diag(S_temp{i});
        if length(sv)<30
            sv = [sv; zeros(30-length(sv),1)];
        end
        s(:,i) = cumsum(sv.^2) ./ sum(sv.^2);
    end
    meanS = mean(s,2);

    subplot(3,1,3)
    plot(meanS,'-','LineWidth',1,'Color',colors{k});%,'Color',colors{k});
    set(gca,'XLim',[0 20],'YLim',[0 1])
    hold on;
    xlabel('Mode Number, n'); ylabel('Cumulative Mode Energy');
    %plot(2,meanS(2),symb(k),'Color',colors{k},'MarkerFaceColor',colors{k});
    
end

legend(actualSpecies{plotInds});

%% Figure S2 plot
plotInds = [2 7 8 9];
ix = 0;
for k=plotInds
    load([specieslist{k} '_' datadate '.mat']);
    ix = ix+1;
    figure(31);
    set(gcf,'pos',[900 100 320 900]);
    subplot(3,1,1);
    dm = size(U_temp{1},1);
    idx = 1:dm:(dm*(length(U_temp)-1)+1);
    Mode1 = cell2mat(U_temp);
    Mode2 = Mode1(:,idx+1);
    Mode1 = Mode1(:,idx);
    std1 = std(abs(Mode1),[],2);
    sOverL = (0:(length(std1)-1))/(length(std1)-1);
    med1 = median(Mode1,2);
    shadedErrorBar(sOverL,med1,std1,'lineProps',{'-','Color',colors{k}});%{'r-o','markerfacecolor','r'});
    hold on;
    ylabel('Curvature Mode 1, u_1 (a.u.)')
    xlim([0 1]); ylim([-0.5 0.5]); set(gca,'YTick',[-0.4 -0.2 0 0.2 0.4]);
    
    subplot(3,1,2);
    std2 = std(abs(Mode2),[],2);
    med2 = median(Mode2,2);
    shadedErrorBar(sOverL,med2,std2,'lineProps',{'-','Color',colors{k}});%{'r-o','markerfacecolor','r'});
    hold on;
    ylabel('Curvature Mode 2, u_2 (a.u.)')
    xlabel('Flagellar Arclength, s/L')
    xlim([0 1]); ylim([-0.5 0.5]); set(gca,'YTick',[-0.4 -0.2 0 0.2 0.4]);
    
    for i = (kfolds(1)*kfolds(2)):-1:1
        sv = diag(S_temp{i});
        if length(sv)<30
            sv = [sv; zeros(30-length(sv),1)];
        end
        s(:,i) = cumsum(sv.^2) ./ sum(sv.^2);
    end
    meanS = mean(s,2);

    subplot(3,1,3)
    plot(meanS,'-','LineWidth',1,'Color',colors{k});%,'Color',colors{k});
    set(gca,'XLim',[0 20],'YLim',[0 1])
    hold on;
    xlabel('Mode Number, n'); ylabel('Cumulative Mode Energy');
    %plot(2,meanS(2),symb(k),'Color',colors{k},'MarkerFaceColor',colors{k}); 
    
end

legend(actualSpecies{plotInds});

%% Figure 2 panels c and d
plotInds = [1 9 4 10 3 6];

ix = 0;
for k=plotInds
    flip1 = strcmp(specieslist{k}, 'Bull');%strcmp(specieslist{k}, 'Bull');% || strcmp(specieslist{k},'Zebrafish')
    flip2 = strcmp(specieslist{k},'Zebrafish');
    switch12 = strcmp(specieslist{k},'Human-low');%strcmp(specieslist{k},'Human-low');
    
    if flip1
       f1 = -1; 
    else
       f1 = 1;
    end
    if flip2
       f2 = -1; 
    else
       f2 = 1;
    end
    
    ix = ix+1;
    load([specieslist{k} '_' datadate '.mat']);
    
    figure(33);
    set(gcf,'pos',[900 100 320 600]);
    subplot(2,1,1);
    dm = size(U_temp{1},1);
    idx = 1:dm:(dm*(length(U_temp)-1)+1);
    Mode1 = cell2mat(U_temp);
    if switch12
        Mode2 = Mode1(:,idx);
        Mode1 = Mode1(:,idx+1);
    else
        Mode2 = Mode1(:,idx+1);
        Mode1 = Mode1(:,idx);
    end
    std1 = std(abs(Mode1),[],2);
    sOverL = (0:(length(std1)-1))/(length(std1)-1);
    med1 = median(Mode1,2);
    shadedErrorBar(sOverL,med1*f1,std1,'lineProps',{'-','Color',colors{k}});%{'r-o','markerfacecolor','r'});
    hold on;
    plot(sOverL(2),med1(2)*f1,symb(k),'Color',colors{k},'MarkerFaceColor',colors{k}); 
    ylabel('Curvature Mode 1, u_1 (a.u.)')
    xlim([0 1]); ylim([-0.7 0.7]);
    
    subplot(2,1,2);
    med2 = median(Mode2,2);
    std2 = std(abs(Mode2),[],2);
    shadedErrorBar(sOverL,med2*f2,std2,'lineProps',{'-','Color',colors{k}});%{'r-o','markerfacecolor','r'});
    hold on;
    plot(sOverL(2),med2(2)*f2,symb(k),'Color',colors{k},'MarkerFaceColor',colors{k}); 
    ylabel('Curvature Mode 2, u_2 (a.u.)')
    xlabel('Flagellar Arclength, s/L')
    xlim([0 1]); ylim([-0.7 0.7]);
    
    for i = (kfolds(1)*kfolds(2)):-1:1
        sv = diag(S_temp{i});
        if length(sv)<30
            sv = [sv; zeros(30-length(sv),1)];
        end
        s(:,i) = cumsum(sv.^2) ./ sum(sv.^2);
    end
    meanS = mean(s,2);

specieslist{k}
median(std1)/(max(med1)-min(med1))
median(std2)/(max(med2)-min(med2))

end

%legend(actualSpecies{plotInds});

